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1. INTRODUCTION 



Source confusion has been a long-standing problem in the astronomical history. In the previous 
formulation of the confusion problem, sources are assumed to be distributed homogeneously on the sky. 
This fundamental assumption is, however, not realistic in many applications. In this work, by making use 
of the point field theory, we derive general analytic formulae for the confusion problems with arbitrary 
distribution and correlation functions. As a typical example, we apply these new formulae to the source 
' confusion of infrared galaxies. We first calculate the confusion statistics for power-law galaxy number 

counts as a test case. When the slope of differential number counts, 7, is steep, the confusion limits 
becomes much brighter and the probability distribution function (PDF) of the fluctuation field is strongly 
' distorted. Then we estimate the PDF and confusion limits based on the realistic number count model 

for infrared galaxies. The gradual flattening of the slope of the source counts makes the clustering effect 
rather mild. Clustering effects result in an increase of the limiting flux density with ~ 10 %. In this 
' case, the peak probability of the PDF decreases up to ~ 15 % and its tail becomes heavier. Though the 

effects are relatively small, they will be strong enough to affect the estimation of galaxy evolution from 
number count or fluctuation statistics. We also comment on future submillimeter observations. 
, Subject headings: methods: statistical — galaxies: statistics — infrared: galaxies — large-scale 

I ' structure of universe — submillimeter — surveys 

cn ■ 
p 

Astronomical source counts often suffer from a problem that multiple sources are located in a single beam of the 
^ ^1 observational instrument. The obtained number counts will then be different from true ones, because the apparent 
I position and flux of the sources are changed by blending of other, usually fainter sources. This is called the source 
O ' confusion, and the resulting measurement error is referred to as the confusion noise. 

^ I In its early stage, confusion problem has been investigated and formulated mainly by radio astronomers in relation 
^ ' to the so-called fluctuation analysis or P{D) analysis (Scheuer 1957; Hewish 1961; Bennett 1962; Murdoch et al. 1973; 
Condon 1974; Condon & Dressel 1978; Wall 1978; Wall et al. 1982). This analysis has been immediately applied and 
developed further in X-ray (e.g., Scheuer 1974; Franceschini 1982; Barcons & Fabian 1990; Barcons 1992; Barcons et al. 
1994) and infrared astronomy (e.g.. Hacking & Houck 1987; Franceschini et al. 1989; Oliver et al. 1997; Matsuhara et al. 
2000; Lagache & Puget 2000; Miville-Deschenes, Lagache, & Puget 2002; Fricdmann & Bouchet 2003). Now confusion 
problem also becomes important in submillimeter cosmology and high-precision astrometric missions (cf. Hogg 2001). 

Cosmologists compare model predictions with observed source counts to extract information of the evolution of galaxies 
and other objects (e.g., Franceschini et al. 1991; Guiderdoni et al. 1998; Takeuchi et al. 1999; Hirashita et al. 1999; Frances- 
chini et al. 2001; Takeuchi et al. 2001a, b). The confusion noise often prevents us from direct comparison between measured 
counts and predictions. Especially, recent advent of infrared (IR) and submillimeter facilities may have stimulated the 
discussion on the confusion problem. Compared with other wavelengths, we have relatively small aperture telescopes in 
the IR, mainly because ground-based observations are almost impossible, and we inevitably need cooled space telescopes. 
In the submillimeter the antenna is relatively large up to ~ 10-30 m, but the long wavelength also makes the single-dish 
surveys confusion-limited by a large diffraction. 

Generally, the measurement error of the source flux blurs the true source counts. This problem was originally pointed 
out by Eddington (1913) and deeply considered in subsequent studies (e.g., Eddington 1940; Bennett 1962; Refsdal 1969; 
Murdoch et al. 1973; Hogg & Turner 1998), but there still remains some unsolved issues, especially when the error is 
dominated by the confusion noise. Recently Monte Carlo simulation has become popular to evaluate confusion effect (e.g., 
Bertin et al. 1997; Eales et al. 2000; Scott et al. 2002), but such an approach is sometimes not easy to interpret, and 
often time consuming. Today, many large survey projects will be performed or completed soon, and a general analytical 
prescription for the confusion is desirable. 

Fundamental assumption for the formulation of the confusion problem is that the sources are distributed homogeneously 
on the sky (e.g., Scheuer 1957; Condon 1974; Franceschini et al. 1989). However, this assumption is obviously not realistic 
in many applications, e.g., stars in the Galaxy, or galaxies in the Universe. A straightforward attempt to take into account 
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the source clustering is to integrate over the spatial correlation functions along the line of sight. This approach has been 
adopted by several works (e.g., Barcons & Fabian 1988; Barcons, Fabian, & Carrera 1998; Burigana & Popa 1998). 
Among them, Toffolatti et al. (1998) presented comprehensive results closely related to the present issue. Since it requires 
the knowledge of three-dimensional correlation functions, while it is a convenient method for theoretical predictions, 
obscrvationally it is quite rare that wc can obtain spatial information at the early phase of a survey. Hence, wc need a 
method to evaluate the confusion only from the projected information of source clustering. In this line, Barcons (1992) 
has pioneered the methodology to tackle the riddle, but his results were restricted to a few simple cases. Since then, little 
analytical progress on this topic has been made up to now. 

In this work, by making extensive use of the methods of the theory of point process, we show the general analytic 
formula for the confusion problems with arbitrary distribution and correlation functions. In Section 2, we first consider 
the fluctuation of unclustered sources, and then present the general formulae for inhomogeneous and clustered source 
distributions. Based on the new formulae, we reconsider how to treat the confusion noise in Section 3. In Section 4 
we formulate higher-order clustering of galaxies, which will be included in the general confusion problems. Wc focus on 
far-infrared and submillimeter galaxies. Section 5 is devoted to our conclusions. Factorial moments and cumulants are 
introduced in Appendix A. In Appendix B we derive the first- and second-order cumulants of clustered field. We present 
a general expression for the rule-of- thumb to avoid confusion in Appendix C. Glossary of the symbols are provided in 
Appendix D. 

2. STATISTICS OF THE FLUCTUATION FIELD 

Statistics of fluctuation field is a fundamental tool for the confusion problem. We first formulate the fiuctuation of 

confusion noise in line with classical works (Scheucr 1957, 1974), and then we translate it to the language of the theory of 
point fields (Daley & Vere- Jones 2003). By this method, we can straightforwardly extend our formulation to the confusion 
problem of clustered sources (cf. Barcons 1992). 

2.1. Fluctuation Generated from Sources without Clustering 
2.1.1. Classical derivation 

We define a probability density function (hereafter PDF) of a source having a flux S' G [S,S + dS], p{S) dS, i.e., 



POO 

/ PiS)dS = l. (1) 



We set h{x) = h{9, (p), the beam pattern (normalized to miity at the beam center), and s{x) = Sh{x), the system response 
from the intrinsic flux S. 

The total signal response we observe, I{x), is described as 

oo 

I{x) = J2 SnKx - Xn) , (2) 

n=0 

where Sn is the flux of the n-th source, and Xn is the angular position of the source. For the following discussion we need 
to deflne a flux signal that consists of exactly N sources, In{x), 

N 

In{x) = ^Snh{x - Xn) ■ (3) 

Consider an ensemble of the number of sources in a unit solid angle, ng, £ = 1,2,..., and let the mean of rii is n. 
Here, 91(5') = np{S) gives the familiar differential number count of the considered sources. Then, from the no clustering 
assumption, the probability of observing exactly N sources in a beam is given by the Poisson process 

P. = (^e--^ (4) 

where fib is the solid angle of the beam. 

Then let us consider the PDF of the value of I{x), f{I),^ 

f{I) = Prob {I' e [1,1 + dl)} 

OO 

= Prob {exactly TV sources lie in a beam} Prob {In{x) G / + dl)} 

N=0 

oo 

= ^PNgN{l), (5) 

N=0 

where gN{I) denotes the PDF of /jv(i)- 

■^This function corresponds to P{D) in radioastronomical terminology. 



Under a certain regularity condition, a PDF is uniquely characterized via its Laplace transform (LT)^. We define the 
LTs of /(/) and gN{I), ^f{t) and Cgf^{t), respectively, 



(6) 
(7) 

(8) 



/OO 
e-"f{l)dl, 
-OO 

/OO 
e-"gi,{l)dl, 
-OO 

where E [ • ] represents the expectation value of a random variable. Then we obtain 

OO 

JV=0 

Hence, concrete expression of Cgj^{t) is our next step to have the functional form of /(/)• 

We consider random variables s„ = Snh{x — x') and their summation over n, In{x) = Ylfn=o Then the PDF of a 
signal / with N summands is represented by the following recursive convolution: 

5jv+i(/)= / gN{I-s')g{s')ds' , (9) 
Jo 

gi{l) = g{l). (10) 

Since the LT of a convolution of two functions is a normal product of the LTs of them, the PDF of s, 

g{s) = Prob {s' G [s, s + ds)} 

gives 



where 



Then we have 



/OO 
e-''g{s)ds. 
-OO 



(11) 
(12) 

(13) 



JV=0 



m 



N 



N 



^ m 



(14) 

The rest wc should consider is the exact form of Cg{t). The profile of a single source is expressed as s{x) = Sh{x). 
Since the LT of a random variable s is a statistical average of e^*", we obtain 



= - / / 

f^b iOb Js 



-tSh(x) 



p{S)dSdx. 



(15) 



Here wc used a simplified symbol = /J^q, i.e., integration with respect to 5 over the range of [0,oo]. We use this 
notation for some other variables in the following. 

Substituting Equation (15) into Equation (14), we finally obtain 



= £-1 



expn 



Jfib \-Js 



-tSh{x) 



p{S)dS-l 



dx 



where C ^ [■] represents the inverse Laplace transform. By the formula for moments 

M.-E[7^]=(-l)'=^ 

and cumulants (reduced moments) 



Kk = (-1) , / 



(16) 
(17) 

(18) 



^If we use it instead of t, we obtain a characteristic function (CF). While the CF is also quite common in physical studies, here we use LT 
to refer literatures in mathematics and statistics. 
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(e.g., Stuart & Ord 1994) we obtain the following beautiful expressions for the cumulants of / (the simplest result of 
Campbell's theorem: Campbell 1909; Rice 1944). 



Ki = E[I] — fj,i — n 



fib JS 



Sh{x)p{S) dSdx = 
K2 = E [(/ -Kif]= [ ( S'^h{xf^{S) dSdx , ... 

jQh Js 



Sh{x)m{S) dSdx , 



and generally, 



S''h{x)''^iS) dSdx . 



(19) 
(20) 

(21) 



2.1.2. Alternative derivation via point field theory 

The position of galaxies, stars, or other point sources is expressed as a point in the field under consideration. A 
mathematical technique to treat such a filed of points is called the theory of point process, and has a long history (e.g., 
Stoyan & Stoyan 1994; Stoyan, Kendall, & Mccke 1995; Daley & Vcrc- Jones 2003). The theory provides us very powerful 
tools for the problems that we consider here. Thus, we translate the above heuristic derivation of /(/) with the language 
of point process theory. Considering the problem along this line enables us a straightforward extension of our discussion 
to the case of clustered sources. 

The fluctuation filed I{x) is expressed as 



I{x) = ^ Snh{x - x„) 



h{x-x')S{x')M{dx') , 



(22) 



where M{A) of Borel sets ^ e is the so-called 'counting measure', which represents the number of points in the set A, 
and S{x') is a fictitious stochastic process that takes a value S{x'^) = Sn at each point x'^ (Daley & Vere- Jones 2003). 
Here we identify the celestial sphere with a real plane R^. We define a random measure A4{A) as 



M{A)= [ I{x)dx 

J A 

= [ [ h{x-x')S{x')N{dx') 

= Sn j h{x - Xn)dx . 



dx 



(23) 



Here we introduce a Laplace functional Lm\X\ with respect to the random measure M., 



exp 



exp 



- / X{x)M[dx) 



X{x)l{x)dx 



(24) 



It is a Laplace counterpart of the characteristic functional of a random field, both of which are often used in particle 
physics, turbulence theory, and structure formation theories in the Universe, etc. (e.g., Vlad et al. 1994; Prisch 1995; 
Szapudi & Szalay 1993; Matsubara 1995)^. We observe 



M2 



X{x)M(dx) 



R2 



X{x) 



dx 



Here, 



exp 



= 2_] / X{x)h{x — Xn)dx . 

/ X{x)h{x — Xn)dx 

JR2 



-Sn / X{x)h{x — Xn)dx 
^In statistical and particle physics, Cj^^ [X] is often expressed in the form 



Z{Xn) . 



(25) 



(26) 



in relation to the path integral formulation (e.g.. Pry 1984a; Szapudi & Szalay 1993; Matsubara 1995). Here ^I] is the assigned probability 
functional for an overall configuration of the field I{x). 
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In Equation (26), >C«r!['] represents the Laplacc-Sticltjes functional with respect to the PDF of 5„, p(>S'„), i.e., 



[ X{x)h{ 



3C CC ) diJC 



/ exp —Sn / X{x)h{x — Xn)dx 

Js L Jr^ 



p{Sn)dSn . 



(27) 



Combining Equations (24), (25), and (26), we have an important relation 



Cm[X] = 1[E 



exp 



-Sn j X{x)h{x — Xn)dx 



= J{Z{Xn) 



= G\Z\. (28) 

The last step is the definition of the probability generating functional (PGFL) of a point field. G [Z] (cf. Balian & Schaeffer 
1989; Daley & Vere-Jones 2003). This shows an important fact: the fluctuation fleld blurred by a beam proflle is expressed 
in terms of the PGFL of the original point field. To be exact, we can describe the observed fluctuation field as an explicit 
functional of the true point field. 

In general, G [Z] has some useful expansions with respect to the factorial moments, factorial cumulants, and other 
related summary statistics of the point field (Vlad et al. 1994; Kerschcr 2001; Daley & Vere-Jones 2003). The most familiar 
statistic for astrophysical studies may be the correlation function (or equivalently, normalized factorial cumulants) . Then, 
G [Z] can be expressed in the following form (for the proof, see, e.g.. Ma 1985): 



\nG[y + l]=f2^ J ■■■ Jy{xi)---y{xk)c^k]dxi---dxk 



= E 

fe=l 



y{xi) ■ ■ ■ y{xk)wk dxi--- dxk , 



(29) 



where C[k] denotes the factorial cumulant of the point field (see Appendix A), and Wk = Wk{xi, ■ ■ ■ ,Xk) is the angular 
fc-point correlation function of the point sources. This relation will be used to extend our formulation to the clustered 
point sources in §2.3. For a Poisson field, by its definition, all the higher order (fc > 2) factorial cumulants vanish, and 
we obtain 



lnG[3; + l] =n / y{x)dx. 



Hence, by substituting Z — y + I in Equation (30) and doing some algebra, we obtain 

G [Z] = e**/B2[-^('=)-il«''= . 

Now returning back to Equation (24), 

jCj^[X] = G[Z] = e"/«2[^W-ild<«= . 
In order to obtain the formula for the local process in a beam area Ob, we set the test function Z{x) as 

Z*{x) = l-[l-Z{x)]Ia^, 

where is an indicator function of a set A, 

'1 X € A, 



x^ A. 



(30) 

(31) 
(32) 
(33) 

(34) 



Substituting Equation (33) into Equation (32) yields 



Recall that 



we obtain 



Z{x„) = £ 



/ X{x)h{x — Xn)dx 



(35) 
(36) 



Cm[X]= e^Y>f^' j {^^<n / X{x')h{x' — x)dx' — 1> dx 

= expn / ( / e-^/"^ x{x')h{x'-x)dx' ^(^g^^g - l] dx . 



(37) 
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-tSh{a 



By setting the test function X{x') = t, we have 

jCf{t) = exp h I 
Lis 

= expn / / e-*^^^''^p{S)dS - 1 



p{S)dS - 1 



dx . 



Inverting the Laplace transform, we again obtain the desired result 



dx . 



(38) 



-tSh{x] 



p{S)dS--l 



dx 



(39) 



/(/) = C ^ expn 

which is the same result with Equation (16). 

2.2. The Case of Inhomogeneous Poisson Point Field 

Sometimes we face a problem that the points distribute locally Poisson but the intensity m is spatially inhomogeneous, 
i.e., depends on the position on the sky. The projected stellar density distribution of the Milky Way might be described 
in this way (e.g., Chandrasekhar & Miinch 1950). 

If the typical angular scale of the variation of h is larger than the typical beam size, we only have to divide the sky into 
patches with the variation scale and derive the PDF of the fluctuation /(/) in each patch by Equation (16). 

On the other hand, the variation scale of n is comparable or smaller than the beam size, we should properly treat the 
variation within a beam. In such case n is expressed as n{x) at the scale of our interest. The intensity in a set A, N{A), 
is obtained by integrating h{x) over A, as 



N{A) 



n{x)dx , 



(40) 



(Cressie 1993, pp. 650-652). In this case the joint probability for some positions are still Poisson, and the independence 
still holds. It leads to the following formula 



Cf{t) = exp /" ( e-'''''''(^)p(5')fiS' - 1 
Jn^, VJ s 



nib i-'s 

The same as the above, the PDF /(/) is obtained by inverting the LT. 



n{x)dx . 



(41) 



2.3. Fluctuation Generated from Clustered Sources 

We now turn to the case of the point source with significant clustering. Clustering of the sources makes the confusion 
effect more severe, because the variance of the source number is larger than that of Poisson fluctuation. The importance 
of clustering in the fluctuation and confusion problems has already pointed out as early as the end of 1960's by Refsdal 
(1969), but he did not provided a quantitative answer in that paper. The fluctuation analysis of sky brightness including 
the effect of clustering was first presented by Barcons (1992) for the study of the X-ray background. The central tool for 
his analysis was characteristic functional of the field I{x). We use its equivalent, Laplace functional here, and derive the 
formula in mathematically rigorous way, in parallel with the discussion in §2.1.2. 

Again we start with the fluctuation field I{x) = X^^g ^nh{x — a;„), but this time x„ are not distributed at random 
on the sky, but have a certain correlation with each other. The statistical properties of the field are characterized by 
Cm['^]- Since the derivation of Equation (28) does not depend on the clustering property of the source point field, we can 
also apply Equations (28) and (29). Since the point field has non- vanishing correlation functions for clustered sources. 
Equation (29) reads 



\nG[Z] 



E 



xi)-l]--- [Z{xk) - l]wk{xi, ■■■ ,Xk)dxi--- dxk ■ 



therefore 



G[Z] 



En 

k=l 



kl 



[2{xi) -!]••• [2{xk) - l]wk{xi, ■■■ ,Xk) dxi ■ ■ ■ dxk ■ 



The same as the unclustered field, we set Z* = 1 — (1 



G[Z] 



e^pY.^J---J[Z{x,)-l]...[. 



Z)lQy^, and obtain 

1] • • • [Z{xk) - l\wk{xi, ■■■ ,Xk)dxi--- dxk ■ 



(42) 



(43) 



(44) 



fib X • ■ ■ X fib 



Substituting Equation (36) yields 



■Ca^ [A"] = exp — / • • • / TT "I £«n / X{x)h{x - Xj)dx - 1> Wk{xi,- ■ ■ ,Xk)dxi ■ ■ ■ dxk ■ 



ObX-'-xOb - 



(45) 
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Again by letting X{x) = t, we have 



CO -I. 



£/(t)=exp^ — 



k=l 



ObX-'-xfib 



oo _u 



Ell 



nhX-'-xOh 



n « 



Wk{xi,--- ,Xk)dxi---dxk 



p{Sj)dSjWk{xi, ■ ■ ■ , Xk) dxi ■ ■ ■ dxk ■ 



Laplace inversion gives the final general formula for the PDF of intensity fluctuation, /(/). 
We obtain cumulants of the PDF of I by the same procedure as Equations (19)-(20): 



Ki 



K2 



Qb JS 



Shix)miS) dSdx , 
S^h{xfm.{S) dSdx 



Qb -'Oi 



2 

n / s,h{^,)misi)ds, 



W2 dXidX2 , 



(46) 



(47) 



(48) 



and so on (sec Appendix B). It is a natural result that ki is the same as that of unchistcrcd sources. On the other hand, 
K2 involves an additional term with respect to W2 = W2{xi,X2), which causes the fluctuation field overdispersed compared 
to that of unclustered sources. 



3. RECONSIDERATION OF THE CONFUSION PROBLEM 

The concept of confusion noise and confusion limit is closely related to the fluctuation fleld /(/). Condon (1974) 
formulated the confusion limit for a power-law source number counts. His formulation was then extended for general 
number counts by Franceschini et al. (1989). Takeuchi et al. (2001b) have made some improvement for precise calculation 
by the method. We should note that the confusion noise affects source counts up to the flux much larger than that of 
confusion limit (Murdoch et al. 1973). We reconsider the confusion problems in this section. 



3.1. Confusion Noise 

First we review the problems related to the confusion noise. Eddington (1913) proposed an important concept that 
existence of any statistical noise makes the observed source counts biased upward. Now this well-known phenomenon is 
called the 'Eddington bias'. He considered the case that the noise is Gaussian. After that, a large amount of efforts has 
been made to treat and correct the bias by subsequent authors (e.g., Eddington 1940; Bennett 1962; Murdoch et al. 1973). 
The outline of these ideas is briefly summarized in Murdoch et al. (1973), in comparison with their own method. 

Problems with Gaussian noise is beautifully formulated by Murdoch et al. (1973). Their setting assumes that the 
dispersion of the error does not depend on the source flux. Even in such a simple case, Murdoch et al. (1973) showed 
that the confusion noise can be very severe for sources with low signal-to-noise ratios. They pointed out that when the 
confusion noise dominate the error, the total noise cannot be Gaussian, and gave a qualitative discussion on the problem. 
They claim to use a Monte Carlo approach to evaluate in that case. 

However, as we will discuss in the subsequent subsections, we have an expression for the fluctuation field, and we can 
utilize the information for evaluating the effect of the confusion noise in general cases. Murdoch et al. (1973) properly 
distinguish between the problems of confusion and blending: the former generally means that the source flux is affected 
by fainter sources, which cannot be detected by the considered instrument individually, and the latter indicates that the 
source flux is affected by fainter sources, which is bright enough to be detected as a source if they were located far from 
the central brightest source. The term 'confusion limit' is related to the former phenomenon, and consequently, confusion 
limit flux is weaker than the noise flux for bright sources. 

How to treat these two effects in a unified manner? Of course they occur in a same way, hence there is no distinction 
in principle. Consider an ideal situation that there are no detector noise, photon noise, read-out noise, nor any other 
cosmological diffuse background radiation with structure. Then the noise intensity distribution is the fluctuation PDF 
itself. 

Here we consider a source with true flux 5*. The situation is schematically described in Figure 1. The source is inevitably 
affected by other sources in a same beam. The total noise intensity in the beam is /. If / < S, then the source is observed 
with flux S + I (see the middle panel of Figure 1). On the other hand, if / is dominated by a bright source S' > S, 
then the source is observed as an additive noise of the brighter source with flux S' > S (the bottom panel of Figure 1). 
However, Murdoch et al. (1973) pointed out that there is a possibility that the signal / consists of two or more sources 
with fluxes fainter than S. The probability of such multiple source blending is obtained by the same manner as discussed 
in §2. Again we should consider sj = Sjh{x — Xj), under the condition 

k 
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Fig. 1. — A schematic description of the confusion noise problem. Top panel: The source with flux S is embedded in the fluctuation field. 
Middle panel: The source is observed with flux 5 + / if the total noise intensity in the beam is /. Bottom panel: the source is observed as 
an additive noise of the brighter source with flux S' > 5 if / is dominated by one or more bright sources S' > S (dashed line). On the other 
hand, if / consists of some sources fainter than S, I should be merely added to S as a noise even if it is larger than S. 



and the probability of having a noise intensity / produced by k sources is expressed by a convolution [Equation (9)]. For 
/ to be a mere noise, there must not be sources brighter than S in a beam. Hence, we should derive the conditional 

probability of having / caused by k sources under the condition Sj < S, j = I, ■ ■ ■ ,k, and then sum up for all ks. It may 
seem to be a complicated problem, but it can be simplified by Laplace transform. By restricting the range of integration 
from to /S in Equation (15) we have the LT of the conditional PDF as 

-I k 

1 



Oh 



(49) 



Then we have the LT of the desired total conditional PDF, /(/; S) 

jC^{t; S) = expn 



as 



dx 



for unclustered point sources, and 

Cf{t; 5) = exp^^ 



fe=i 



k 

n 



e-ts'M''j)p{S'j)dS'j 



1 



Wk dxi ■ ■ ■ dxk 



(50) 



(51) 



for clustered sources. 

Hence, in order to obtain the confused number counts from the true ones, we should convolve the fluctuation distribution 

/(/) itself at / < S', and /(/; 5) at / > S. This consideration naturally explains the fact that the confusion noise flux 
is stronger than the confusion limit flux. Since we have a variety of background in the image data in a real situation, 
we fit and subtract from the obtained flux, the noise is not necessarily positive. Thus in practice, the dependence of the 
confusion noise distribution on flux S causes only a weak variation along S under the existence of other kind of noise, 
because of the convolution with other noise. 
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3.2. Confusion limit 

We next formulate the relation between the beam size dh and the source confusion limit''. The basic procedure is 
to estimate the limit signal strength as an upper cutoff of the integration in the second-order cumulant formulae in 
Equation (20) or Equation (48) by iterations. Here we express the second-order cumulant formulae as a function of signal 
strength s: 

' S'^h{xfmiS) dSdx 



K2 



= I S 

for unclustered sources, and 

S'^h{xf^{S)dSdx+ Iff 
iob Jn^, Jsi 



h{x) 
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for clustered sources. We utilize these formulae for deriving confusion limits in the following. 

3.2.1. Confusion Limit for power-law number counts 
We begin our discussion with the case that the number count is described by a power-law, according to Condon (1974)^: 



m{s) = aS-'' . 

For unclustered sources, we obtain the confusion limit flux to a cutoff signal Sc from Equation (52) as 



(54) 
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where QeS is the so-called effective beam size, defined as 
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Taking the square root of Equation (55) and setting Sc = aa as often used, we have 



If the beam pattern is described by a Gaussian 



3-7 



(aOeff) 



1/(7-1) 



HO) 



where is the FWHM of the beam, we get Condon's analytic formula for the a-a confusion limit. 
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It is convenient to relate these expressions and the empirical rule of thumb for the confusion limit. We discuss the relation 
between the above formulae and the empirical rule of thumb in Appendix C. 

Now wc turn to the same problem for clustered sources. For clustered sources, the confusion limit can be obtained in 
similar way by substituting Equation (54) as follows: 



Jo JQb 



h{x) 

= r s2 / ah{ey-^ 

Jo JQi, 



dx 
h{x) 



ds + 



siS2 [ [ "yi 

Jflb JOb 



Si 



s-T eded(j) + 



eff 



3-7 



4-' + {^is) 



SlS2 



,2(2-7) 



h{xi) 
[ [ a''[h{ei 

Jcib v/Qb 



S2 



h{X2) 



W2 



dxi dx2 
h{xi) h{x2) 



dsids2 



\j,/iaM7-i/ ^-7 0id9id(j)i 92de2d(j}2 , , 

)h{02)\ (S1S2) W2 — TTTTl TTTTl — dSidS2 



h{0l) h{e2) 



(60) 
(61) 



where we defined the following quantity 

{^ls)= I I [h{ei)h{92)r~^w2eideid(f>ie2de2d(t)2 

In general, we cannot solve Equation (60) analytically, hence we should calculate it numerically. 

''Hereinafter we assume that h{x) = {6,(t>) = h{6) for simplicity. This is not an essential assumption for the subsequent discussions. 

*Here 7 is the power-law exponent of differential number counts. This is the same convention with Pranceschini (1982). Note that Barcons 
(1992) uses the same character 7 as the exponent of cumulative number counts, and consequently the related expressions are apparently different 
from those in Barcons (1992) in terms of 7. Hogg (2001) also use the cumulative count slope for his numerical study of the confusion errors. 
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3.2.2. Confusion Limit for general number counts 

Pranceschini et al. (1989) presented an iterative formula for the confusion limit of the general number counts. Here we 
see the result for a Gaussian beam [Equation (58)]: 
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where rj = e'*i"2(e/eb)^ g^^^j^ ^^le upper cutoff that corresponds to the beam area. We again set Sc = aa, then we have 

(63) 
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Thus we obtain a general relation between the beam size and the confusion limit as 
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Same as the case for power-law number counts, when source clustering takes place, we cannot have a simple formula, 
but we can still extract some information from the expression. Again by setting the upper cutoff Sc in the integrals in 
Equation (53), we have 
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■ ■ ■ ,Xk), similar numerical 



Hx) Jo 

If the beam profile is a Gaussian, and some specific functional form is obtained for = wi.{xi, 
computation can be performed in parallel with Equation (62). 

4. APPLICATION: GALAXY CLUSTERING AND CONFUSION 
4.1. Hierarchical Ansatz 

As we have discussed above, in order to treat the chistering of the sources, a set of prescribed angular correlation 
functions to some reasonable order arc required. Here we focus on the galaxy clustering as a typical example of the 
related issues. Of course, we can handle any point sources as far as we have some knowledge of their clustering. 

Property of higher-order galaxy clustering still remains a matter of debate (for a comprehensive overview, see Peebles 
1980), and various models have been advocated in order to describe the correlation function of galaxies. Generally, a 
hierarchy of the correlation functions continues to infinite order, hence we should introduce a certain closure relation to 
cut the sequence. 

The hierarchical Ansatz is often employed for the correlation functions, both from phenomenological basis (e.g., Balian 
& Schaeffer 1989) and theoretical considerations of BBGKY equations (e.g., Davis & Peebles 1977; Fry 1984b; Yano & 
Gouda 1997). The Ansatz claims that the fc-point correlation is described by a product of — 1 two-point correlation 
functions. These studies consider the three-dimensional correlations, but we can relate them to two-dimensional ones in 
a simple way (Peebles 1980; Szapudi & Szalay 1993). We apply the following hierarchical relations for angular galaxy 
correlation on the sky (Meiksin, Szapudi, & Szalay 1992): 

(66) 



Wk{xi,--- ,Xk) =qk 51 W2{Xri,Xr2)---W2{x. 



rk-2' Xrk-i) ) 



where Qk is a numerical factor to determine the strength of clustering, and the summation J2rj trees taken over all tree 
topologies of graphs between points (Fig. 2). For the details of the expansion, see Fry (1984a, b). 



4.2. Correlation of Infrared Galaxies 

Now we focus on the distribution of infrared galaxies, to have some insight to the confusion in the forthcoming infrared 
and submillimeter galaxy surveys. Today, many large infrared and submillimcter missions are planned or in progress, 
and their source confusion limits for galaxies have been calculated for each facility (e.g., Ishii, Takeuchi, & Sohn 2002; 
Dole, Lagache, & Puget 2003). They all assumes the random distribution of galaxies on the sky. In order to obtain the 
confusion limits more precisely, we must take into account the clustering of infrared galaxies properly. 

By using the density moment technique, Meiksin, Szapudi, & Szalay (1992) estimated the coefficients qs—qs in Equa- 
tion (66) for IRAS 1.2 Jy sample. We use these values up to A; = 4 to approximate the clustering: = 1.25 and 54 = 2.19. 
The two-point angular correlation function is 



w{xi,X2) = W {\xi - X2I) = W{6i2) = 
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Fig. 2. — The graph representation of the hierarchical expansion of the correlation functions. 



and (5 = 0.79, = 0?36 for the IRAS 1.2 Jy sample. We can obtain the clustering of any deeper surveys via scaling 
relations of the correlation functions with the characteristic depth of the survey, (Peebles 1980): 

W2{ei2) = d:^w^idJi2) , (68) 

W3{ei2, 023, 03l) = d-^wl{d,e^2, dj23, dJ3l) , (69) 
W4,{6i2, 623, ^34, Oil, Ol3, 624) = ^U>4((i,(?i2, d^923, d*6'34, d*^41, d^.613, d^924) , (70) 

where superscript represents that it is evaluated at 5 = 1.2 Jy at A = 60 /xm, i.e., W2{9) = (^^/0?36)~"■^^, and is 
the relative characteristic depth of the survey, defined as = (1.2 [Jy]/S' [Jy])^/^. We observe that the contribution of 
higher-order clustering rapidly decreases with increasing depth of the survey. This dilution of clustering is more effective 
for higher-order correlations. 

Now we obtain the two-point correlation function W2 {0) at an arbitrary flux 5 through Equation (68) as 

W2{0) = d-^'+^^w°2{9) = d:'-''wU0) 
'1.2[Jy]N-^-^^/^ 



(71) 



Consequently, from Equations (69) and (70), we have 
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and 
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Fk;. 3. — The power-law number count model as a model of infrared galcixy counts. These model counts are the integrated ones, hence 
the power-law indices are represented by 7 — 1 = 1.5, 2.0, and 2.5. 



4.3. Results 

4.3.1. Power-law number counts 

In order to evaluate how strong the correlation structure affects the confusion-related quantities, we first calculate these 
statistics based on power-law number count models [Equation (54)], 

91(5') = aS-'' , 

with the power-law indices 7 — 2.5, 3.0, and 3.5, i.e., the indices of the integrated counts are 7 — 1 = 1.5, 2.0, and 2.5, 
respectively. Here, consider a survey at A = 60 iim. The power-law number counts are normalized so that they have the 
same counts with the IRAS QMW galaxy survey at a flux S = 0.9 Jy (Rowan-Robinson et al. 1991). It is observationally 
known that the bright end of the counts is well approximated by Euclidean, because cosmological and evolutionary effects 
are both negligible. Actually, IRAS QMW galaxy counts shows a good fit to the Euclidean slope (slope index 7 = 2.5 in 
differential counts), and we can safely fix the slope of the model 7 to be 2.5 above the flux S > 0.9 Jy. We estimated the 
confusion-related statistics for the cases with and without clustering. We assume a telescope with an aperture of 70 cm 
with an ideal Airy function as a PSF, and the detection limit of the instrument is 50 mJy. 

As discussed above, the clustering of galaxies depends on the detected flux: brighter sources have stronger clustering 
on the sky, and the clustering gradually becomes weaker toward fainter flux. This makes the exact formulation for 
the statistical characteristic of the two-dimensional galaxy distribution prohibitively difficult (see Barcons 1992), and 
unfortunately, an exact mathematical theory to treat this problem has not fully established yet. In order to calculate 
the confusion limit with clustering, we approximate the clustering of the whole sample galaxies evaluated at a 'fiducial' 
flux, S'fid, instead of the flux-dependent clustering in gradual way. Fluctuation consists not only of the detected sources 
but also of unresolvable sources fainter than detection limit, in principle, toward infinitesimally faint flux. Therefore, 
the fiducial flux is fainter than the point source detection limit. We assumed that the fiducial flux is proportional to 
detection limit. Based on this assumption, we calibrated the fiducial flux empirically so as to reproduce the IRAS 
confusion limit (cr ~ 20 mJy: Hacking & Houck 1987; Lonsdale & Hacking 1989; Berlin et al. 1997). We foimd that the 
relation Sua = 5iini/40 can be used to evaluate the average clustering strength of the IRAS galaxy sample. Namely, all 
the galaxies with S > 5fid are approximated to have the same correlation functions [w2{0) = d^^W2{dfidO), etc.] and 
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The 5-a confusion limits for infhared calaxjes with i'ower-law number counts. 
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Fig. 4. — The probability density function (PDF) of the fluctuation /(/) for a power-law number count model with a slope index 7 = 1.5. 
We present the PDF in linear scale (Left panel) and in logarithmic scale (Right panel) . The dotted lines show the fluctuation distribution for 
unclustered (randomly distributed) sources, whereas the solid lines represent that for the sources with clustering. 
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Fig. 5. — The PDF of the fluctuation for a power-law number count model with a slope index 7 = 2.0. In contrast to the case of 7 = 1.5, 
the effect of clustering is drastic. 



ignore the contribution from the sources with S < S'fid- Of course, exact treatment of this part remains to be theoreticaUy 
improved. 

We calculated the confusion Umits from the power-law model counts. Without clustering, the confusion limits would be 
3.1 mJy, 22.9 m.Jy, and 138 mJy for the indices 7 — 1 = 1.5, 2.0, and 2.5, respectively. If we take into account the effect of 
clustering properly, they become 3.8 mJy, 120 mJy, and 9.33 Jy, for 7 — 1 = 1.5, 2.0, and 2.5, respectively (see Table 1). 
Thus, the confusion-limit flux get larger, and the effect of clustering strongly depends on the count slope. Especially, if 
the cumulative count slope 7 — 1 exceeds 2.0, its effect will be catastrophic. 

We also obtained the PDF of the fluctuation intensity for infrared galaxies under the same assumptions. Figure 4 shows 
the PDF with the integral number count slope index 7 — 1 = 1.5. The peak probability decreases by clustering. The 
clustering also result in a broadening of the PDF (Barcons 1992). It also should noted that the probability of finding a 
very low intensity at a certain position (closely related to 'the void probability') increases when clustering takes place. 

The effect of clustering becomes very strong when the slope index approaches 7 — 1 = 2.0, as presented in Figure 5. In 
this case, we have already found that the clustering makes the confusion limit much shallower than that of unclustered case. 
It means that the rms of the confusion noise is large. Actually, we found that the broadening of the PDF by clustering is 
very strong, and result in a large variance of the fluctuation PDF (Fig. 5). In summary, the effect of clustering can be very 
severe if the number counts are described by a single power-law. The above result shows that the relative contribution of 
clustering to the total fluctuation compared with Poisson component becomes larger for large 7. 

These arc consistent with the results reported by Toffolatti ct al. (1998) for dusty galaxies at IR wavelengths. As 
mentioned above, their approach is more model-oriented than the present work, based on the simple model of the spatial 
two-point correlation function of galaxies. Their central aim was to estimate the power spectrum of the fluctuation in 
the submillimeter and radio background, and they considered only the two-point correlation, i.e., second-order statistics. 
Since confusion limit is the second-order quantity, hence both methods can be used, and yield consistent estimates, though 
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our method requires only projected information. 

Our fornndation is fundamentally based on the PDF, and it is the most suitable method of calculating it from the 
measured projected information (see §2). In order to extend Toffolatti et aL's method to calculate the PDF, it is necessary 
to model the bispectrum, trispectrum, and higher-order spectrum or their equivalents (for thorough formulation, see 
Matsubara 2003). It is a theoretically interesting but challenging problem, and again we see that their method is suitable 
to theoretical models. 

4.3.2. Realistic number counts for infrared galaxies 

The observed number counts of infrared galaxies cannot be described by a power-law function. Deep galaxy surveys 
by ISO revealed some striking features of the infrared galaxy counts (see, e.g., Pranceschini et al. 2001; Pearson 2001; 
Takcuchi et al. 2001a; Lagache, Dole & Puget 2003). Galaxies detected at 90 /im or 170 show a rapid increase in 
number counts at 5 ^ 50-100 mJy, with cumulative slope 7 — 1 > 2.5. There have never been any model that succeeded 
in reproducing the extremely steep number-count slope of ISO sources perfectly, even though many elaborate models have 
been presented to date. 

But such a rapid increase of infrared galaxy cannot be continued toward a fainter flux level, because the integrated 

infrared flux from galaxies must not exceed the observational constraint of the cosmic infrared background (for the details, 
see, Takeuchi et al. 2001a; Pranceschini et al. 2001; Hauser & Dwek 2001). This makes a change of the slope of the observed 
number counts of galaxies. 

We show some calculations of the confusion limits for such realistic counts based on an empirical model of Takeuchi et al. 
(2001a). The parameters are taken from Takeuchi et al. (2002) and Ishii, Takeuchi, & Sohn (2002). The calculation method 
is the same for the case of power-law number counts, except that we should use Equation (65) instead of Equation (60). 

We obtained the confusion limits of 1.8 mJy for unclustered sources and 2.0 mJy for clustered sources for a model number 
count at A = 60 /xm. They are 11.4 mJy and 12.3 mJy at 90 /zm, and 105 mJy and 120 mJy at 170 /xm. We summarize 
these limits in Table 2. Thus in the far-infrared, source clustering does not affect the confusion limit significantly. This 
is explained as follows: Since the cosmic expansion makes the slope flatter toward fainter fluxes, the count slope becomes 
sub-Euclidean at the flux level we consider here, in spite of the strong galaxy evolution. The effect of clustering is 
dominated by the faintest clustered sources, hence the decrease of differential counts makes the effect rather small. 

Then, let us consider the PDF of the fluctuation for realistic number counts. We calculated the PDF for wavelengths 
A = 60,90, and 170 /xm. Figures 6, 7, and 8 show the PDFs of the fluctuation in these wavelengths. Just as seen in the 
power-law count models, the peak probability decreases up to ~ 15 %, and the tail of the PDF becomes heavier. We also 
find an increase of the probability of finding a small flux density. The clustering effect may not seem very severe, because 
present data generally suffer from more serious noise sources and/or systematics. However, we must keep this effect in 
mind when we estimate galaxy evolution from the fluctuation analysis of the data that will be obtained by ASTRO-F and 
SIRTF. 

Remind that we assumed the constant clustering strength in this analysis. One may expect a smaller clustering at 
higher redshift, but biased galaxy formation scenario predicts even stronger clustering for high-redshift dusty galaxies 
(e.g., Ishii, Takeuchi, & Enoki 2003). Thus, we should take care of the effect in future infrared data analysis. 

4.3.3. Comments on submillimeter galaxies 

In the submillimeter wavelengths, confusion effects are very serious in single-dish observations (Blain, Ivison, & Small 
1998; Takeuchi et al. 2001b). Although some large interferometric facilities are under construction, the confusion problem 
will remain an important issue that should be properly considered for future airborne or space submillimetric observational 
projects such as BLAST (Balloon-borne Large-Aperture Sub-millimeter Telescope)^ and Herschel Space Observatory^^ . 
As seen above, we need a two-point correlation function of dusty galaxies, whose reliable observational estimate has not 
been available yet (Scott et al. 2002; Almaini et al. 2003; Borys et al. 2003). For the estimation of the galaxy confusion 
effect in future submillimeter large surveys (see, e.g., Takeuchi et al. 2001b, and references therein), accurate estimation 
of the clustering properties is desired. The measured confusion limit from SCUBA survey data was found to be shallower 

"URL: http: //chilel .physics . upenn.edu/blastpublic/index.shtml. 
^'^URL: http: //astro . estec . esa.nl/First/ . 

Table 2 

The 5-a confusion limits based on the model number counts of infrared galaxies. 



Wavelength 5-a confusion limits [mjyj 



[/xm] 


Random 


Clusterinj 


60 


1.8 


2.0 


90 


11.4 


12.3 


170 


105 


120 



15 




Fig. 6. — The PDF of the fluctuation /(/) for a realistic number count model at 60 //m based on Takeuchi et al. (2001a). We present the 
PDF in linear scale (Left panel) and logarithmic scale (Right panel). The dotted lines show the PDF for unclustered (randomly distributed) 
sources, whereas the solid lines represent the PDF for clustered sources. 




Fig. 7. — Same as Figure 6, except that it is a PDF at 90 //m. 




Fig. 8. — Same as Figure 6, except that it is a PDF at 170 fim. 



than that expected from classical formula (D. H. Hughes 2003, private communication), hence the effect of clustering 
may take place. If we assume that the SCUBA galaxies have the same clustering property with IRAS galaxies, clustering 
increases the confusion limit from 20 % to 100 %, depending on the adopted assumption for the spectral energy distribution 
of galaxies. 

Most of the submillimeter galaxies are dusty vigorous starbursts (e.g., Pranceschini et al. 2001; Takeuchi et al. 2001a, b; 
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Totani & Takcuchi 2002), and the dominant constituent of the cosmic submiUimeter background radiation. Hence, the 
statistical structure of the submiUimeter background carries crucial information of the spatial distribution of high-redshift 
starburst galaxies (e.g., Peacock et al. 2000). Our analysis plays an important role as an interface between observed 
submiUimeter data and sophisticated cosmological models (e.g., Magliocchetti et al. 2001; Perrotta et al. 2003). 



The source confusion is a long-standing problem of the astronomical history, going back to the works of Eddington. 
Fundamental assumption of the formulation is that the sources are distributed homogeneously on the sky. However, this 
assumption is not realistic in many applications, e.g., stars in the Galaxy, or galaxies in the Universe. Clustering increases 
the confusion effect to some extent, but it is not an easy task to formulate the effect in the confusion problem. By using 
a model for spatial correlation functions, Toffolatti et al. (1998) showed a useful result for the contribution of clustering 
to the fluctuation of the background at submiUimeter and radio wavelengths. Observationally, however, the method to 
evaluate the confusion only from two-dimensional information of clustering is desirable. In this line, only the work of 
Barcons (1992) has made attempt to tackle the difficult problem and given solutions to a few simple cases. 

In this work, by making extensive use of the methods related to the point field theory, we derived general analytic 
formulae for the confusion problems with arbitrary distribution and correlation functions. With these formulae, we 
can obtain the statistical properties of the confusion noise caused by the sources with some prescribed two-dimensional 
correlation structure. 

Based on the general formulae, we first calculated the confusion limits from the power-law galaxy number counts as 
a test case for the analysis of infrared galaxies. We considered an infrared facility with an aperture of 70 cm and an 
ideal Airy PSF. For galaxy clustering, we adopted the hierarchical Ansatz to calculate higher-order correlation functions 
used in our formulation. Without clustering, the confusion limits would be S.lmJy, 23 mJy, and 138 mJy for the indices 
7 — 1 = 1.5, 2.0, and 2.5, respectively. If we take into account the effect of clustering properly, they become 3.8 mJy, 
120 mJy, and 9.9 Jy, for 7 — 1 = 1.5, 2.0, and 2.5, respectively. We also obtained the probability density function 
(PDF) of the fluctuation intensity, /(/). We found a broadening of the PDF by clustering, i.e., increase of the variance, 
corresponding to the increase of the confusion limit flux density. When 7 approaches to 2.0, the clustering effect becomes 
very severe for a power-law number counts. 

Then we estimated the PDF and confusion limits based on the realistic number count model of Takeuchi et al. (2001a), 
under the same hypothetical infrared facility. At 60 fim, we obtained confusion limits of 2.0 mJy for infrared galaxies. At 
90 ^m and 170 /xm, the confusion limits are 12.3 mJy and 120 mJy, respectively. These estimates are not much different 
from the limits calculated without clustering [l.BmJy (60 /im), 11.4mJy (90 /im), and 105 mJy (170 /im)]. The gradual 
flattening of the number count slope toward the fainter flux densities, which is observed for infrared galaxies, makes the 
clustering effect very small, even though the slope is very steep at 5 ~ 10~^-10~^ Jy. As for the PDF at these wavelength 
bands, the peak probability also decreases up to ^ 15 %, and the tail of the PDF becomes heavier. The clustering effect 
may not seem very severe, because present data generally suffer from more serious noise sources and/or systematics. 
However, we must keep this effect in mind when we estimate galaxy evolution from the fluctuation analysis of the data 
by ASTRO-F and SIRTF. 

We will also apply our method to future submiUimeter large surveys. Most of the submiUimeter galaxies are dusty 

vigorous starbursts, and the dominant constituent of the cosmic submiUimeter background radiation. Our analysis plays 
an important role as an interface between observed submiUimeter data and sophisticated complex cosmological models. 
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much. We offer our gratitude to Dave Hughes, Motohiro Enoki, Taihei Yano, Haruhiko Ueda, Seiji Ueda, Takuji Tsujimoto, 
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We present the deflnition and some formulae for the important statistical quantities, factorial moments and factorial 
cumulants here. First we consider a probability generating function (PGF), G{u), 



where A'': a nonnegative integer-valued random variable. Here we define the factorial power. For any integers n and k, 



5. CONCLUSIONS 



APPENDIX 



A. FACTORIAL MOMENTS AND CUMULANTS 
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Factorial moments of are defined as 
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where 

p„ = Prob{A^ = n} . 

Then the k-th factorial moments, m^i.] are obtained as the coefficients of the following Taylor expansion, 



G(l + ^) = 1 + ^ 



■m[k]V 
k\ 



(A4) 
(A5) 



Note that wc assumed the convergence of the series. 
Factorial cumulants, c^k], are similarly calculated by 



\nG{l + v) = l + Y, 



k\ 



The first few relations are as follows: 



C[i] = = E [TV] , 

C[2] = "^[2] - m^i] , 

C[3] = TO[3] - 3m[2]m[i] + 2mfi] 



(A6) 

(A7) 
(A8) 
(A9) 

(AlO) 



For a Poisson distribution, pjv = n^e "/A?'!, 

C[i] = n, C[fe] =0 (fc > 2) . 

These summary statistics are extensively used in the field of structure formation, especially for counts-in-ccUs analysis 
(e.g. Balian & Schaeffer 1989; Szapudi & Szalay 1993; Szapudi et al. 1995, 1999). For some more details, see e.g., Vlad et 
al. (1994) and Kerscher (2001). 

B. DERIVATION OF FIRST- AND SECOND-ORDER CUMULANTS FOR GENERAL FIELD 

Here we present how to derive the first- and second-order cumulants for the fiuctuation field of clustered sources. We 
start with the LT [Equation (46)]: 



£/(t) =exp^ 



OO _ 



fe=i 



k\ 



Ob X ■ • • X fib 

Then, by Equation (18), we observe 

dln£f{t) 



k 

n 



-tSjhixj) _ ^ 



p{Sj)dSjWk{xi,- ■■ ,Xk) dxi ■ ■ ■ dxk ■ 



Ki = (-1) 



E 

k=l 



dt 



k\ 



t=0 

d_ 

dt 



Qb X ■ • • X Ob 



n 



Here we define 



U{xj;t)= I 



-tSjh{Xj) _ 



^-tSjh{xj) _ ^ 



^{Sj)dSj > Wk dx\ ■ ■ ■ dxk 



m{Sj)dSj. 



t=0 



then Equation (Bl) becomes 

OO ^ 



k=l 



k\ 



E 



flbX-'-xOh 



dU{xf,t) 
Jt 



U{xi;t)---U{xr,t)---U{xk;t) 

^ V ' 

(k — 1) terms 



Wk dxi ■ ■ ■ dxk 



where • means that the term is omitted in the product. We see that 

U{xj; 0) = / (1-1) ^{Sj)dSj = , 

JSi 



(Bl) 



(B2) 



(B3) 



(B4) 



therefore, among the terms summed over k = 1, • • • , oo, all those with k > 2 vanish and only k = 1 term remains. Thus 
we obtain 



Jn 



dU{x;t) 



fib dt 



dx 



-J J 

;0 ^Hb JS 



-Sh{x)e-*^''<^''A m{S) dSdx 



t=o 



Ob -fs 



Sh{x)m{S) dSdx . 



(B5) 
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Next, we consider the second-order cumulant, H2'- 



t=o 



a \ - 1 

fe=i 



dt 



Ob X • ■ ■ X fib 



dt^ k\ j "' 
''-^ nbx---xab 



^ A;! 



/s=l 



Ed 
Jt 



fibX-'-xOb 



1=1 



n 

dU{xe; t) 
Jt 

dU(xe;t) 
Jt 



^-tSjh{Xj) _ 



^{Sj)dSj > Wk dxi ■ ■ ■ dxk 



U{xi;t) ■ ■ ■U{xr,t) ■ ■ ■U{xk;t) 



j=i-fe 



(/c— 1) terms 

Wfc dXi • • • dXk 



Wk dxi ■ ■ ■ dXk 



t=0 



t=0 



= E 



fe=i 



fibX-'-xOb 



fr dt/(a;,;^) dU{x„,;t) ^^^^ 
2^ 11 U{xr,t)+ — [[ U{x,;t) 



Wk dxi ■ ■ ■ dxk 



df^ ■>•■»■" \ , dt dt 

Among the first term in the square brackets [ • ] of Equation (B6) , all the summands with k > 2 vanish by substituting 
t = 0. On the other hand, the second term in [ • ] exists only if fc > 2. Further, all the summands with k >3 vanish by 
setting t = 0. Considering all the above items together, we obtain 



K2 = 



Ob 



Jn 



d^U{x-t) 

d'^U{x\t) 
d^ 



dx 



dx 



fib "'fi 



dU{xi;t) dU{x2]t) dU{x2]t) dU{xi;t) 



Hb 



fib -IS 



4/ / 



dt dt 
dU{xi;t) dU{x2]t) 



dt 



dt 



dt dt 

W2{Xi,X2)dXidX2 



W2{xi,X2)dXidX2 



t=0 



t=0 



+ 



Ob "'SI 



52/i(a;)2e-*^'*(=^)0T(5) dSdx 

Sih{xi)e-*^'^^'''^m{Si) dSi [ S2h{x2)e-*^^''^''^^m{S2) dS2 

1 •'S2 



W2{xi,X2)dXidX2 



Qb JS 



S^h{xfm{S) dSdx 



fib -T! 



Si 



Sih{xi)m{Si) dSi / S2h{x2)m{S2) dS2 



W2{xi,X2)dxidx2 ■ (B7) 



C. THE EMPIRICAL RULE OF THUMB FOR THE CONFUSION LIMIT 



The '1/30 sources per beam' rule is well known to observational astronomers as an empirical convention to avoid 
confusion effects. This problem was first theoretically addressed by Franceschini (1982). Here we extend his discussion 
and present a numerical result for the case of Airy beam. According to the notation of Hogg (2001), we use the notation 
(s/b) as the number of sources per beam: 



(s/6)(S') = l)b / m{S')dS' . (CI) 
Js 

If we would like to detect a-a sources safely, where a is the confusion limit [Equation (57)], the rule of thumb is expressed 



as 



{s/b)c = flb ^{S)dS, 

J aa 



(C2) 

where the subscript C represents that (s/b) at the confusion limit. For power-law source counts, we can obtain a useful 
analytic expression for {s/b)c as follows: 



{s/b)c = Ob / aS-'^dS 

J aa 



7^ 



(aa) 



1-7 



an 



7-1 

-a 



a 



3-j X (l-7)/(7-l) 



7-1 
fib / 3 " 7 



3-7 
3^ 



(afieff) 



(l-7)/(7-l) 



(af^eff)"' 



(C3) 
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Fig. 9. — The rule of thumb for the confusion limit of 5-iT sources. Dotted line represents the number of sources per beam (s/6)c for the 
confusion limit of a Gaussian beam, as a function of the slope 7 of the differential source counts. Open squares are (s/6)c for an Airy beam. 
These estimates are calculated for unclustered sources. 

This corresponds to the result obtained for the general case by Franceschini (1982) in his Equation (14). For a Gaussian 
beam, Equation (C3) becomes simpler. If we set Ob = ''"(ecc)^, where ctg = ^b/2\/2 ln2 denotes the standard deviation 
of the Gaussian beam profile, we have 

Here e is a numerical factor that depends on the details of the source detection algorithm and beam profile. We regard 
e = 2 as a reasonable choice. Hogg (2001) defined fl^ to be ttctg^, which is different from our convention by a factor of e. 
Then, we observe the following formula 

(^/&)c= ''^^~^^ (C5) 

This formula clearly shows the dependence of the confusion rule of thumb, {s/b)c- the steeper the slope is, the more 
severe the confusion becomes. This trend has already known to observational astronomers empirically, and confirmed by 
the comprehensive work of Hogg (2001). We tabulate the rule of thumb when we take a = 5 and e = 2 in Table 3. 

Generally, an ideal beam pattern of an optical instrument is often described by an Airy function. We calculated the 
rule of thumb for an Airy beam by numerical integration. This is shown in Figure 9. The numbers of sources per beam 
for Gaussian and Airy beam show an excellent agreement at 7 > 1.7. However, we find a slight deviation from Gaussian 
value at 1.3 < 7 < 1.7 for (s/6)c of an Airy beam. This may stem from the structure of the Airy function. 

D. GLOSSARY OF MATHEMATICAL SYMBOLS 

We tabulate the mathematical symbols used in the main text in Table 4, for the readers' convenience. We tabulate the 
equations or sections by which these symbols are defined or explained. 



Table 3 

The rule of thumb for the confusion limit of 5-a sources, when the beam is Gaussian. 



Index 

7 


Sources per beam 

i^/bh 


0.5 


0.20 


1.0 


0.16 


1.5 


0.12 


2.0 


0.08 


2.5 


0.04 


2.9 


0.008 
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Table 4 

Glossary of mathematical symbols. 



Symbol Definition Equation 



Prob {Event} Probability that a certain event occurs (5) 

E [X] Expectation value of a random variable X (6) 

Ia Indicator function of a set A (34) 

R2 Real plane §2.1.2 

Q]^ Beam area of an instrument on the sky (33) 

Q^ff Effective beam size, i.e., beam area weighted by h(dy'^^ (55) 

(^eff) Squared beam area weighted by h(Oi)'^^^ , h{92)'^^^, and W2{0i2) (60) 

a Proportional constant for a power-law source number counts (54) 

P Slope index of the power-law angular correlation function of galaxies (67) 

7 Exponent of a power-law source immber counts (54) 

6,cf),0i,92, ■ ■ ■ Angular coordinate or separation on the sky (55) 

Oij Angular separation between objects i and j (67) 

0b FWHM of the beam pattern (58) 

8(1 Correlation scale length of angular correlation 67 

Kfe fc-th order cumulant of a random variable (18) 

fc-th order moment of a random variable (17) 

cr Confusion limit flux to a cutoff of a-a level (57) 

cr(sc) Confusion limit flux to a cutoff signal Sc (55) 

(TG Standard deviation of a Gaussian beam profile Appendix C 

G[Z] Probability generating functional (PGFL) of a point field (32) 

/(as) Intensity of a signal at a position x (2) 

In{x) Intensity of a signal at a position x, consisting of exactly N sources (3) 

Cf{t) Laplace transform (LT) of f{t) (6) 

C j{t\ S) LT of / under the condition that any St < S for any i = 0, 1, 2, • • • (50) 

Cg{t) LT of g{t) (13) 

Cg^it) LToigN(t) (7) 

'^aki^'t ^) °^ 9k under the condition that any Si < S for any i = 1, ■ ■ ■ ,k (49) 

[u] Inverse Laplace transform of a certain function u{t) (16) 

Cm l'^] Laplace functional of a test function X with respect to a measure M (24) 

£511 [U] Laplace-Stieltjes functional of U with respect to p{S) (26) 

M{A) Random measure a random field with respect to an area A (23) 

N Measured number of sources in a beam (3)-(5) 

Af{A) Counting measure of a point field in an area A (22) 

01(5) Differential source number counts with a flux density S §2.1.1 

S Detected flux density of a source §2.1.1 

S„ Detected flux density of a source with label n §2.1.1 

Slim Detection limit flux for point sources §4.3.1 

5fid Fiducial flux to approximate the average correlation of a sample §4.3.1 

X Test function, with which a value is fixed for a functional (24) 

y Test function of a functional (29) 

Z{xn) Laplace-Stieltjes functional of J^2 X{x)h{x — Xn)dx (26), (28) 

Z*{xn) Test function Z(a3„) but locally restricted within a beam (33) 

a Flux density level measured by standard deviation of the noise §3.2.1 

C[fe] k-th factorial cumulant (29), (A6) 

d* Relative characteristic depth of a survey with respect to S = 1.2 Jy (68)-(70) 

/(/) Probability density function (PDF) of the signal intensity of the image (5) 

g{I) PDF of Ii, i.e., signal produced by a single source (10), (11) 

gN{I) PDF of /jvi i-e., signal produced by exactly A'' sources (5), (9) 

h{x) Beam pattern of the instrument §2.1.1 

m[(i.] fc-th factorial moment (A5) 

n Mean of a source number density on the sky, n£ {t = 1,2, . . .) §2.1.1 

nl'^l fe-th factorial power of a certain integer n (A2) 

Pl^ Probability of obtaining a number A' for a Poisson random variable (4) 

p(5) Probability of finding a source with a flux S §2.1.1 

g/j Coefficient of the hierarchical model for the fc-point correlation function (66) 

s{x) Intensity of the signal at position x, s{x) = Sh{x) §2.1.1 

Sc Cutoff flux level used for calculating confusion limit (55) 

Sn{x) Intensity of the signal at x produced by a single source at Xn §2.1.1 

{■^/b)c Number of sources per beam at the confusion limit (C3) 

{s/b){S) Number of sources per beam to a flux limit S (CI) 

Wkixi, ■ ■ ■ ,xi^) Angular fe-point correlation function (29), (67) 

w^{xi, ■ ■ ■ , Xh) Angular fe-point correlation function evaluated at S = 1.2 Jy (68)-(70) 

X Position where the intensity of the field is measured §2.1.1 

Xn Position of a source with label n §2.1.1 
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